/*******************************************************************************
																			
	DESCRIPTION: 	This code plots the distribution of the estimated parameters
					in the cyclicality regressions.

*******************************************************************************/

clear all
global id_code 114_4
set seed 2110

* B1: Import data from R
use "${data}/114_3_CyclicalityBetas_Full_2006.dta", clear

/***********************************************************************
* B3: Histograms of coefficients from the panel regression
************************************************************************/
* Calculate e^beta_0
gen exp_beta_0_shrunk_ind = exp(beta_0_shrunk_ind)

* Calculate 20th and 80th percentiles
sum exp_beta_0_shrunk_ind
local mean_0: di %9.3f r(mean)
local sd_0: di %9.3f `=`r(sd)''
local var_0: di %9.3f `=`r(sd)'^2'

_pctile exp_beta_0_shrunk_ind, p(20 80)
local p20_0: di %9.3f `r(r1)'
local p80_0: di %9.3f `r(r2)'

sum beta_u_shrunk_ind
local mean_u: di %9.3f r(mean)
local sd_u: di %9.3f `=`r(sd)''
local var_u: di %9.3f `=`r(sd)'^2'

_pctile beta_u_shrunk_ind, p(20 80)
local p20_u: di %9.3f `r(r1)'
local p80_u: di %9.3f `r(r2)'

_pctile beta_t_shrunk_ind, p(20 80)
local p20_t: di %9.3f `r(r1)'
local p80_t: di %9.3f `r(r2)'

* Bootstrap the shrunken coefficients:
bootstrap r(sd), reps(500) seed(2110) verbose: sum beta_u_shrunk_ind

* Plot histograms:
histogram exp_beta_0_shrunk_ind, frequency  width(0.02) fcolor(ebblue*0.5) lcolor(ebblue) ///
	legend(on order( ///
		- "Mean = `mean_0'" ///
		- "Std. Dev. = `sd_0'") ///
		symxsize(*0.5) size(small) cols(1) pos(10) ring(0) justification(right) placement(right) ) ///
	graphregion(color(white)) /// 	
	plotregion(margin(b=0 y=0)) ///
	xaxis(1 2)					///
	ytitle("Frequency", size(9pt)) ///
	xtitle("exp({&beta}{sub:0})", axis(1) size(9pt)) ///
	xtitle("", axis(2))		///
	ylabel(0(2000)8000, format(%11.0fc) labsize(10pt) angle(0)) yscale(titlegap(2)) ///
	xscale(titlegap(2) axis(1)) xlabel(0(0.25)1, axis(1))	///
	xscale(noline axis(2)) ///
	xlabel(none, value axis(2) noticks) ///
	name(dist_beta0, replace)
graph export "${output}/${id_code}_Cyclicality_Distribution_Beta0Shrunken.pdf", as(pdf) replace
		
histogram beta_u_shrunk_ind if inrange(beta_u_shrunk_ind, -0.8, 0.2), frequency width(0.02) fcolor(ebblue*0.5) lcolor(ebblue) ///
	legend(on order( ///
		- "Mean = `mean_u'" ///
		- "Std. Dev. = `sd_u'") ///
		symxsize(*0.5) size(small) cols(1) pos(2) ring(0) justification(right) placement(right) )	 ///
	graphregion(color(white)) /// 	
	plotregion(margin(b=0 y=0)) ///
	xaxis(1 2)					///
	ytitle("Frequency", size(9pt)) ///
	xtitle("{&beta}{sub:U}", size(9pt)) ///
	xtitle("", axis(2))		///
	ylabel(0(2500)17500, format(%11.0fc) labsize(10pt) angle(0)) yscale(titlegap(2)) ///
	xscale(titlegap(2) axis(1))	xlabel(-.8(0.2).2, axis(1))		///
	xscale(noline axis(2)) ///
	xlabel(none, value axis(2) noticks) ///
	name(dist_beta1, replace)
graph export "${output}/${id_code}_Cyclicality_Distribution_BetaUShrunken.pdf", as(pdf) replace				

histogram beta_t_shrunk_ind, frequency fcolor(ebblue*0.5) lcolor(ebblue) ///
	legend(on order( ///
		- "Mean = `mean_u'" ///
		- "Std. Dev. = `sd_u'") ///
		symxsize(*0.5) size(small) cols(1) pos(2) ring(0) justification(right) placement(right) )	 ///
	graphregion(color(white)) /// 	
	plotregion(margin(b=0 y=0)) ///
	xaxis(1 2)					///
	ytitle("Frequency", size(9pt)) ///
	xtitle("{&beta}{sub:T}", size(9pt)) ///
	xtitle("", axis(2))		///
	yscale(titlegap(2)) ///
	xscale(titlegap(2) axis(1))	///
	xscale(noline axis(2)) ///
	xlabel(none, value axis(2) noticks) ///
	name(dist_betat, replace)
graph export "${output}/${id_code}_Cyclicality_Distribution_BetaTShrunken.pdf", as(pdf) replace				

* Compate beta_t and beta_u:
heatplot beta_t_shrunk_ind beta_u_shrunk_ind, ///
	xtitle("{&beta}{sub:U}") ///
	ytitle("{&beta}{sub:T}") ///
	title("Joint Distribution of {&beta}{sub:U} and {&beta}{sub:T} (both shrunken) ") ///
	legend(subtitle("% of sample", size(small))) keylabels(, format(%5.2f)) ///
	graphregion(color(white)) name(heat, replace)
	
graph export "${output}/${id_code}_Cyclicality_Distribution_TvsU_heatplot.pdf", as(pdf) replace				

